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Abstract 

A method is proposed and studied whereby the 
system identification cycle consisting of experiment 
design and data analysis can be repeatedly 
implemented aboard a test aircraft in real time. This 
adaptive in-flight system identification scheme has 
many advantages, including increased flight test 
efficiency, adaptability to dynamic characteristics that 
are imperfectly known a priori , in-flight improvement 
of data quality through iterative input design, and 
immediate feedback of the quality of flight test results. 
The technique uses equation error in the frequency 
domain with a recursive Fourier transform for the real 
time data analysis, and simple design methods 
employing square wave input forms to design the test 
inputs in flight. Simulation examples are used to 
demonstrate that the technique produces increasingly 
accurate model parameter estimates resulting from 
sequentially designed and implemented flight test 
maneuvers. The method has reasonable computational 
requirements, and could be implemented aboard an 
aircraft in real time. 
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Nomenclature 

body axis translational accelerations, g 
expectation operator 
linear system matrices 
acceleration due to gravity, ft/sec 1 2 
altitude, ft 
cost function 
aerodynamic lift force 
Mach number 

total number of sample times 
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p,q,r body axes angular velocities, rad/sec 

Re real part 

R measurement noise covariance matrix 

Tr trace 

u control vector 

V airspeed, ft/sec 

x state vector 

y t output vector at time i At 

z t measured output vector at time i At 

Z z body axis aerodynamic force 

0 zero vector 

a angle of attack, rad 

Sij Kronecker delta 

At sampling interval, sec 

(O angular frequency, rad/sec 

0 /^-dimensional parameter vector 


superscripts 

T 

t 


-1 


transpose 

complex conjugate transpose 
discrete Fourier transform 
estimate 
matrix inverse 


subscripts 

o 

w 


trim or initial value 
wind axes 


Introduction 

Flight testing to collect data for dynamic modeling 
begins with some a priori information about the 
aircraft dynamics. The a priori information typically 
consists of the rigid body equations of motion in 
conjunction with aerodynamic and engine data from 
ground tests. Flight test maneuvers are designed to 
excite the dynamic response of the aircraft, based on 
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the a priori information to produce data for dynamic 
modeling. The maneuvers are then scheduled for 
flight test, executed by either the pilot or an on-board 
computer system in flight, and the data is recorded for 
post-flight analysis. Because of differences between 
flight test results and predictions based on ground 
tests, the maneuver design based on a priori 
information may turn out to be deficient in some way. 
For example, if the input amplitude is too low, the data 
will have low information content, which leads to 
inaccurate model parameter estimates. A similar data 
information deficiency can result if the input 
frequencies are poorly chosen. On the other hand, if 
the input amplitude is too high, the aircraft motion may 
stray too far from the test condition or excite 
nonlinearities, which can invalidate modeling 
assumptions made for both the maneuver design and 
the data analysis. Typically, these problems are 
discovered during post-flight data analysis, resulting 
in a need for additional flight tests. Such additional 
experimentation requires more flight time, more 
engineering time, and more money, which are seldom 
readily available. A block diagram of the entire 
process appears as Figure 1. 

One iteration through the identification cycle 
shown in Figure 1 can easily take months or even 
years in practice. The main difficulty is the large time 
delay that occurs due to the analysis and scheduling of 
additional maneuvers, and the inevitable conflict of the 
requirement for more flight time to collect dynamic 
modeling data with other objectives of the flight test 
program. 

This work examines a different approach to 
obtaining the requisite data for dynamic modeling 
purposes. Figure 2 illustrates the idea. Analysis is 
now done in real time aboard the aircraft, so that 
modeling results from an initial maneuver could be 
used to design the subsequent test maneuver, and so 
on throughout the flight test until acceptable results 
(e.g., specific model parameter accuracies) are 
achieved. The result is an adaptive experiment design 
and data analysis method that is carried out in flight. 
With this approach, the quality of modeling results 
based on the measured flight data would be known 
before the aircraft lands, and no further flights would 
be required for this purpose. At the very least, the 
procedure could provide a high level of confidence that 
sufficient high quality data was collected for the 
modeling task, and could do so in one flight. Flight 
testing in this way would therefore greatly increase 
efficiency and effectiveness, since the complete 
identification cycle would be carried out in flight in an 


2 


automated fashion. In addition, detailed a priori 
information would not be required to design the 
maneuver in this case, since the maneuver design 
starts by using a very simple low amplitude input, and 
then evolves the input to more complex forms with 
appropriate amplitudes based on analysis of data from 
the preceding maneuvers. 

There are many candidate methods for the 
in-flight data analysis and maneuver design, but the 
main requirement that they be simple enough to be 
implemented in real time aboard the aircraft narrows 
the field. In particular, any methods that iterate 
through the recorded data must be eliminated. One 
possibility is to use an extended Kalman filter to 
estimate the model parameters. This approach has 
been described in the literature 1 , along with several 
specific applications to aircraft parameter estimation 
problems 2-5 . There are some problems with this 
approach, however, some of which are related to the 
fact that the constant model parameters are treated as 
additional states, which can therefore vary with time. 
In addition, numerical and convergence problems 
related to the requisite linearization and noise variance 
estimation can be encountered. 

In this work, a simple recursive computation of 
the Fourier transform is used to implement equation 
error in the frequency domain for in-flight model 
parameter estimation. The in-flight maneuver design 
is done by creating simple square wave input forms 
with progressively increasing complexity, based on 
the latest results from the in-flight model parameter 
estimation. Together, these techniques adaptively 
carry out the complete system identification cycle in 
flight, using simple methods that can be implemented 
on modem flight computers. 

The next section gives the problem statement and 
outlines the necessary theory. Following this, the 
in-flight system identification method is applied to a 
simulation example, where a linear truth model is used 
with outputs corrupted by noise similar to that 
observed in flight. The application is identifying an 
accurate model for the rigid body short period 
dynamics of a conventional fighter. Finally, a 
nonlinear simulation is used to demonstrate the 
in-flight system identification procedure. 
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Theoretical Development 

Airplane dynamics can be described by the 
following linear model equations: 


x(t) = Ax(t) + Bu(t) 

(1) 

II 

o 

H 

(2) 

y(t) - Cx(t) + Du(t) 

( 3 ) 

ii 

+ 

ii 

JO 

;N (4) 


The discrete measurement noise vector v t is assumed 
Gaussian with 

E{y t ) = 0 and 2?{v,v/} = *5 /y (5) 

Matrices C, and D in Eqs. (1) and (3) 
contain stability and control derivatives, which are the 
constant model parameters to be estimated from flight 
test data. The input quantities are control surface 
deflections, with output quantities from air data 
( V, a, [i ), body axis angular velocities ( p,q,r ), 
Euler angles ( (j), 0, y / ) , and translational accelerations 
( a x , a y , a z J. Longitudinal and lateral cases can be 

treated separately, with the linear model structure 
shown above resulting from the usual small 
perturbation assumptions 6 . 

The linear model structure given here will be used 
for the in-flight maneuver design and the in-flight data 
analysis. Using the linear model structure keeps the 
in-flight calculations simple, and is adequate for the 
purpose at hand, which is to adaptively design flight 
test maneuvers that produce good flight test data. If 
necessary, more sophisticated modeling and data 
analysis techniques can be applied to the measured 
data post- flight. 

Equation Error in the Frequency Domain 

The finite Fourier transform of a signal x(t) is 
defined by 

x(co) = x(t)e~ J(0t dt (6) 


which can be approximated by 
N - 1 

x((o) » At ^x i e~ j(0ti (7) 

i= 0 

Subscript i indicates the variable value at time iAt, 
and At is the sampling interval. The summation in 
Eq. (7) is defined as the discrete Fourier transform, 


N - 1 

X(oj) = Y, x i e~ JC0ti 

( 8 ) 

i = 0 


x(co) ~ X(co) At 

( 9 ) 


Some fairly straightforward corrections 7 can be 
made to Eq. (9) to remove the inaccuracy resulting 
from the fact that Eq. (9) is a simple Euler 
approximation to the finite Fourier transform of 
Eq. (6). If the sampling rate is much higher than the 
frequencies of interest (a)) , then the corrections are 
small and can be safely ignored. 

Applying the Fourier transform to Eqs. (1) and 
(3) gives 

jo)x(co) = Ax{(o) + Bu(co) (10) 

y(co) = Cx(co) + Du(co) (11) 

When the states, outputs, and inputs are measured, 
individual state or output equations from vector 
Eqs. (10) or (1 1) can be used in an equation error 
formulation to estimate the stability and control 
derivatives contained in matrices A, B,C, and D . 
For the Ath state equation of vector Eq. (10), the cost 
function is 

1 m ^ 

Z n = 1 

( 12 ) 

where A k and B k are the Mi rows of matrices A and 
B , respectively, and x k (a) n ) is the Ath element of 
vector x((D n ). There are m terms in the summation, 
corresponding to m frequencies of interest, and each 
transformed variable depends on frequency. Similar 
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cost expressions can be written for individual output 
equations from vector Eq. (11). Denoting the vector 
of unknown model parameters in A k and B k by 6 , the 
problem can be formulated as a standard least squares 
regression problem with complex data, 



( 19 ) 


and p is the number of elements in parameter vector 6 . 


where 


y e =x c e + £ c 


M **(®i) 

jO) 2 Xk((0 2 ) 


Y c = 


j(O m X k ((O m ) 


(13) 


(14) 




X r (© 1 ) U r («!) 

x t (co 2 ) u t (co 2 ) 

\x T {(0 m ) u T (co m ) J 


(15) 


and £ c represents the equation error in the frequency 
domain. The least squares cost function is 


J= l -{Y c -X c d) f {Y c -X c d) 


(16) 


which is identical to the cost in Eq. (12). The 
parameter vector estimate that minimizes this cost 
function is computed from 8 


d = [Re(xJX c )\ l Re(xjY c ) (17) 


The estimated parameter covariance matrix is 

cov(§) = e{(§ - efe - e) r } = o 2 [Re(x;x c )] 1 

(18) 

where the equation error variance a 2 can be estimated 
from the residuals, 


Recursive Fourier Transform 

For a given frequency, the discrete Fourier 
transform in Eq. (8) at sample time i is related to the 
discrete Fourier transform at time i — 1 by 

X i {(0) = X^yim) + Xj e~i aiht (20) 

where 

e -j(OiAt - e -j(oAt e -j(o(i-\)At (21) 


The quantity e~ JG>At is constant for a given 
frequency and constant sampling interval. It follows 
that the discrete Fourier transform can be computed for 
a given frequency at each time step using one addition 
in Eq. (20) and two multiplications - one in Eq. (21) 
using the stored constant e~ JQ)At for frequency CQ, and 
one in Eq. (20). There is no need to store the time 
domain data in memory when computing the discrete 
Fourier transform in this way, because each sampled 
data point is processed immediately. Time domain 
data from all preceding maneuvers can be used in each 
subsequent analysis by simply continuing the 
recursive calculation of the Fourier transform. More 
data from more maneuvers improves the quality of the 
data in the frequency domain without increasing 
memory requirements to store it. In addition, the 
Fourier transform is available at any time i A t. The 
approximation to the finite Fourier transform is 
completed using Eq. (9). 

Rigid body dynamics of piloted aircraft lie in the 
rather narrow frequency band of 0-1 Hz. It is 
therefore possible to select closely spaced fixed 
frequencies for the Fourier transform and the 
subsequent data analysis. In this work, frequency 
spacing of 0.02 Hz was found to be adequate, which 
gives 50 frequencies evenly distributed on the interval 
[0.02-1.0] Hz for each transformed time domain 
signal. Zero frequency is excluded to remove trim 
values and measurement biases. The number of time 
domain signals to be transformed is usually low (7 or 
less - more if there are many control surfaces), so that 
this approach requires a small amount of computer 
memory that is independent of the time length of the 
flight test maneuvers. 
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Test Maneuver Design 

The specification of the flight test maneuvers 
(equivalently, the flight test inputs) has a major impact 
on the quality of the measured data for modeling 
purposes. The goal is to design an experiment which 
produces data from which model parameters can be 
accurately estimated. This translates into exciting the 
dynamic response modes such that the sensitivities of 
the model outputs to the parameters are high and 
correlations among the parameters are low. Designing 
an experiment which meets these objectives requires 
rich excitation of the system, which is frequently at 
odds with various practical constraints, such as the 
requirement that output amplitude excursions stay 
within specified limits in order to assure the validity of 
an assumed model structure. 

Previous work 9 ’ 10 has shown that inputs 
optimized for parameter estimation experiments can be 
effectively designed using the V -optimality criterion, 
wherein the input design cost J T equals the sum of 
squares of estimated parameter variances: 


Ji = TA 


E\id-d)[d-d^ 


( 22 ) 


The analytical connection between the input time 
history and the achievable accuracy for the model 
parameter estimates based on the measured data is 
detailed in Ref. [9]. For the present purposes, it 
suffices to say that the connection is strong, which 
means that the input time history has a significant 
impact on the accuracy of the parameter estimates 
computed from measured data. The choice of input 
implicitly includes the length of the maneuver. 

Doublet, 2-1-1, and 3-2-1-1 input forms have the 
advantages of easy implementation in flight and simple 
design based on current estimates of modal 
frequencies and steady state (dc) gain. Figures 3, 4, 
and 5 illustrate each input form. 

Flight test results from Ref. [9] demonstrated that 
at low angle of attack, a simple 3-2-1-1 input form is 
roughly 25% less effective than a globally optimal 
square wave input which minimized the V -optimality 
criterion given above. 

The 3-2-1-1 has a very simple design procedure, 
which is: 

1 . Match the frequency of the 2 pulse to the current 
estimate of the natural frequency for the 
dominant oscillatory mode. Note that a single 
pulse represents one-half the period. 
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2. Scale the 3 and 1 pulse widths in proportion to 
the 2 pulse. 

3 . Set the amplitude of the pulses so that output 
amplitudes do not exceed values that would 
invalidate the assumed model structure using the 
current model estimate. For linear dynamical 
systems, this is a simple scaling operation. 

The design procedure is similar for the 2-1-1 
input: 

1 . Select the pulse width so that the frequencies of 
the 2 and 1 pulses bracket the frequency of the 
current estimate of the natural frequency for the 
dominant oscillatory mode. In the current work, 
the 2 and 1 pulses were determined as 4/3 and 
2/3 times the pulse width corresponding to the 
current natural frequency estimate. 

2 . Set the amplitude of the pulses so that output 
amplitudes do not exceed values that would 
invalidate the assumed model structure using the 
current model estimate. For linear dynamical 
systems, this is a simple scaling operation. 

Doublet inputs are sometimes designed to match 
the frequency of the current estimate of the natural 
frequency for the dominant oscillatory mode, but are 
also sometimes designed to have very thin pulse width 
to approximate a two-sided impulse, which 
theoretically contains all frequencies. 

There are other input forms and input design 
methods 10 , but square wave input forms have been 
shown to be simple and effective 9 , so these input 
forms were selected for the flight test input design. 

The test maneuver sequence begins with a simple 
unit amplitude doublet with a one second pulse width, 
which corresponds to a frequency of 0.5 Hz. After 
executing this initial maneuver, the data are analyzed in 
real time using equation error in the frequency domain, 
as described above. The requisite Fourier transforms 
for the parameter estimation in Eq. (17) are computed 
continuously in real time using Eqs. (21), (20), and 
(9). Based on the results from this analysis, a slightly 
more complex input, the 2-1-1, is designed and 
implemented. Input amplitudes and frequencies are 
calculated based on results from data produced by the 
previous doublet maneuver. Analysis of data from all 
preceding maneuvers (i.e., the doublet and 2-1-1) is 
used as the basis for designing amplitude and 
frequency for the slightly more complicated 3-2-1 -1 
input. Re-designs of the 3-2-1 -1 input can continue 
based on the most recent data analysis until specific 
objectives (such as specified accuracy on any or all 
estimated parameters) are met. The sequential increase 
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in the input complexity is designed to correspond with 
the increasing accuracy of the model parameters as the 
maneuvers are executed. Total time for each maneuver 
includes 2-3 seconds of zero input before each test 
input to allow responses to settle to trim values in 
preparation for each maneuver. The delay is 
randomized uniformly on the interval [2,3] seconds to 
enrich the input spectrum. This delay also allows time 
for the data analysis and input design for the next 
maneuver. 


Example 

For longitudinal aircraft short period dynamics, 
the state vector x, input vector u, and output vector y 
in Eqs. (1) and (3) are defined by 

x = [a q] T w = [<5 e ] 

r i t ^ ^ 

y = [a q a z \ 


System matrices containing the model parameters are: 

Z* 


A = 


Za 

1 

N 

M a 

£ 

i 


B = 




(24) 



i 

0 


0 

c = 

0 

1 

D = 

0 


sT 

N 3 

1 


I 

N 

1 


(25) 


The above model assumes a ~ q, so that a effects 
can be subsumed into the Z' and M q derivatives. 
Parameter Z' includes the inertial term, i.e., 

Z' = 1 + Z q . The model also assumes a small trim 
angle of attack, so that Z~-L and a 7 ~a 7 . 

In the first example, perturbation elevator inputs 
were applied to a known linear model to produce 
simulated state and output responses. The simulated 
aircraft is a conventional F-16 6 with forward c.g. 
position (0.2 c ) in straight and level flight at 
10,000 ft, trim angle of attack 7 deg, and Mach 0.37. 
The simulated states and outputs were corrupted with 
20% gaussian random white noise. This made the 
signal-to-noise ratio 5-to-l for each simulated state and 


output measurement. The elevator input was assumed 
to be measured without noise, which is a close 
approximation to reality. The true values of the model 
parameters used to generate the simulated test data are 
given in column 2 of Table 1. Parameter estimation 
was done using equation error in the frequency 
domain applied to the two state equations, as described 
above, with the Fourier transform computed 
recursively. 

The initial doublet input is shown in Figure 3. 
Parameter estimation results based on simulated data 
from this input are given in column 3 of Table 1. In 
general, the parameter estimates are approximate, and 
the standard errors (in parentheses below each 
parameter estimate) are relatively large. Figure 4 is the 

2- 1-1 input design based on the results from the 
doublet maneuver data. The abscissa in Figure 4 
shows that the 2-1-1 maneuver was executed after the 
doublet maneuver. Input amplitude was scaled to 
produce ±2.5 deg of a change from trim, which is 
typical for linear model validity. The input amplitude 
scaling was done by keeping track of the maximum 
absolute a deviation from trim during the doublet 
maneuver, a max (in deg), then scaling the doublet 
input amplitude by the ratio 2.5/ a max . Pulse widths 
were designed according to the procedure outlined 
above. Parameter estimation results in column 4 of 
Table 1 show that the parameter estimates are 
approaching the true values with standard errors 
decreasing. A similar trend is shown in column 5 of 
Table 1, which contains results based on the data 
including the 3-2-1-1 maneuver shown in Figure 5. 
The 3-2-1-1 maneuver was designed based on analysis 
of data from the preceding doublet and 2-1-1 
maneuvers. The abscissa in Figure 5 shows that the 

3- 2-1 -1 maneuver was executed after the 2-1-1 
maneuver. 

Figure 6 shows the trend for the Z a parameter 
with each identification cycle of maneuver design and 
data analysis. The estimated parameter values 
approach the true value with a decreasing standard 
error. Similar plots could be made for the other model 
parameters. 

The technique was then applied to the full 
nonlinear F-16 simulation 6 , with similar results. 

Figure 7 shows the elevator input design for a 
simulated flight test of the F-16 with forward c.g. 
position ( 0.2 c ) in straight and level flight at 
20,000 ft, trim angle of attack 20 deg, and 
Mach 0.27. Figure 8 shows the angle of attack and 
pitch rate response. The same output amplitude 
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constraint of ±2.5 deg in a was used and was 
successfully implemented via the input amplitude 
scaling. 

Parameter estimation results for this case are 
given in Table 2 using the same format as Table 1 . 
The "true" linear model parameter values shown in 
column 2 of Table 2 were obtained from the nonlinear 
simulation using central finite differences and 1% 
perturbations on the state and control variables. 
Performance of the algorithm was similar to that seen 
for the linear case at a lower trim angle of attack. 

Input amplitudes and pulse widths were adjusted 
adaptively and automatically to match the dynamics of 
the aircraft at this flight condition. Figure 9 shows the 
power spectrum of the input. The adaptive input 
design resulted in an input with a broad energy band 
centered at the "true" natural frequency, even though 
the algorithm had no information about the natural 
frequency of the system a priori . 

The nonlinear simulation, the on-line maneuver 
design, and the on-line data analysis were all 
programmed in Matlab and ran on a Sun 200 MHz 
HyperSparc (serial processor, running SunOS 4.1.4) 
roughly twice as fast as real time. 


Concluding Remarks 

The approach to aircraft system identification 
proposed here changes the philosophy of aircraft 
dynamic modeling experimentation from designing test 
maneuvers based on a priori predictions of the 
dynamic characteristics and evaluating the data quality 
post- flight, to an in-flight adaptive approach that relies 
solely on measured flight data from the dynamical 
system to be modeled. The developed method 
implements the identification cycle of experiment 
design, parameter estimation, re-design of the 
experiment based on the estimation results, and so on 
repeatedly, until desired accuracy measures for the 
model parameters are met. Dynamic effects that are 
impossible to predict on the ground before the flight 
could be accounted for in real time by the automatic 
design of the test maneuvers in flight. A priori input 
design is avoided altogether. The procedure has 
reasonable computational requirements and could be 
implemented in flight in real time. 

This work used simulation examples to 
demonstrate the in-flight system identification scheme. 
The method could be used for dimensional or 
non-dimensional parameter estimation, and could also 
be used with general nonlinear models, as long as the 
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model is linear in the parameters. All states and inputs 
must be measured, but this should not be a problem in 
the flight test environment for which the method is 
intended. A more sophisticated algorithm would be 
required for multiple input maneuver design, but 
repeated applications of the present method to 
individual inputs, one at a time, might prove to be 
adequate. 

The on-line data analysis in the frequency domain 
has the advantage of automatically removing trim 
values and measurement biases from the data because 
zero frequency is omitted from the Fourier 
transformation. In addition, the data is automatically 
filtered because only the specific frequencies 
corresponding to the dynamic motion of interest are 
included in the Fourier transformation. 

Ultimately, the algorithm could be packaged as a 
subroutine to be included in flight control computer 
software. Since the algorithm is adaptive and requires 
no a priori analysis for the input design, it could be 
called on to execute appropriate test maneuvers for any 
flight condition throughout the flight envelope. This 
capability has obvious attraction for flight envelope 
expansion, flight control system design validation, 
aerodynamic parameter estimation, and simulator 
updates. 
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Table 1 Linear Simulation Results, 
a 0 = l deg, h 0 = 10,000 ft, M 0 = 0.37 


Maneuver 

Parameter 

True 

Value 

Doublet 

2-1-1 

3-2-1-1 

z a 

-0.600 

-0.791 

(0.054) 

-0.629 

(0.030) 

-0.599 

(0.022) 

7' 

0.950 

1.106 

(0.046) 

0.987 

(0.021) 

0.961 

(0.016) 

Zo 

°e 

-0.002 

0.0010 

(0.0016) 

-0.0015 

(0.0008) 

-0.0016 

(0.0006) 

K. 

-4.300 

-3.598 

(0.137) 

-4.118 

(0.094) 

-4.193 

(0.043) 

M q 

-1.200 

-1.662 

(0.116) 

-1.178 

(0.067) 

-1.184 

(0.030) 

M 5e 

-0.090 

-0.104 

(0.004) 

-0.088 

(0.002) 

-0.089 

(0.001) 


Table 2 Nonlinear Simulation Results. 

a Q 

= 20 deg, h 

0= 20,000 ft, M 0 = 0.27 

Maneuver 

Parameter 

True 

Value 

Doublet 

2-1-1 

3-2-1-1 


-0.260 

-0.611 

(0.140) 

-0.266 

(0.040) 

-0.274 

(0.029) 

7' 

0.954 

1.103 

(0.101) 

0.886 

(0.036) 

0.911 

(0.028) 

Zc 

°e 

-0.0007 

0.0023 

(0.0021) 

-0.0014 

(0.0009) 

-0.0009 

(0.0007) 

M a 

-1.761 

-1.613 

(0.108) 

-2.004 

(0.039) 

-1.878 

(0.033) 

M q 

-0.611 

-0.769 

(0.078) 

-0.591 

(0.035) 

-0.614 

(0.032) 

M Se 

-0.0341 

-0.0381 

(0.0016) 

-0.0333 

(0.0009) 

-0.0342 

(0.0008) 


8 


American Institute of Aeronautics and Astronautics 




Figure 1 Conventional System Identification Cycle Figure 2 In-Flight System Identification 
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Figure 8 Nonlinear Simulation Response 



Figure 9 Input Power Spectrum 
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